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Abstract 

In this paper, we study in detail, both analyticahy and numerically, the dynamical 
properties of the triangle map, a piecewise parabolic automorphism of the two- 
dimensional torus, for different values of the two independent parameters defining 
the map. The dynamics is studied numerically by means of two different symbolic 
encoding schemes, both relying on the fact that it maps polygons to polygons: in 
the first scheme we consider dynamically generated partitions made out of suitable 
sets of disjoint polygons, in the second we consider the standard binary partition of 
the torus induced by the discontinuity set. These encoding schemes are studied in 
detail and shown to be compatible, although not equivalent. The ergodic properties 
of the triangle map are then investigated in terms of the Markov transition matrices 
associated to the above schemes and furthermore compared to the spectral proper- 
ties of the Koopman operator in L^(T^). Finally, a stochastic version of the triangle 
map is introduced and studied. A simple heuristic analysis of the latter yields the 
correct statistical and scaling behaviours of the correlation functions of the original 
map. 
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1 Introduction 



The research of conservative dynamical systems has a long and fruitful history. 
From the discovery of the exponential sensitivity on initial conditions by Poincare 
the main part of research directed its attention to such chaotic situations. The 
Kolmogorov- Sinai entropy is one of the basic invariants of dynamical systems, which 
divides them into two classes: (i) the chaotic or hyperbolic systems, with positive 
dynamical entropy, and (ii) non-hyperbolic systems, with zero dynamical entropy. 
A particular subclass of non-hyperbolic systems has been discovered recently, ex- 
hibiting properties common mostly to chaotic systems such as diffusion, ergodicity 
and mixing. The study of such systems is important in order to understand the 
origin of stochastic-like-dynamics in non-hyperbolic dynamical systems - and in 
particular in non-equilibrium statistical mechanics - and the effects on their quan- 
tum counterparts. A phy sically relevant example in this class is the triangle map 
( Casati and Prosen . 20001 ). As we will see below, this is a two parameter (a, (3) fam- 
ily of two dimensional maps on the torus. They are related to the returning map 
of a point particle moving inside an elongated rectangular triangle billiard, namely 
for the correspondence to be accurate one of the angles of the triangular billiard 
has to be small. Roughly speaking a is related to one angle in the triangle and 
j3 is introduced in order to generalize the map and to make the dynamics richer. 
For the detailed constructio n of the map and its relat i on to the triangular billiards 
we refer to original articles (ICasati and Prosenl . Il999l . l2000f). For a more complete 
description of the dynamics in polygonal billiards see e.g. (iMasur and Tabachnikovl . 
2OO2I : iTabachikovl . Il9 95) . Some k nown generic p roperties of dynamics in polygonally 
shaped billiards are reviewed in (lGutkinl . [T996h . The triangle map : — i> is a 
one-to-one transformation (automorphism) of a two-dimensional torus = [0, 1)^ 
written in additive notation (mod 1) as 



P) = {q + P + a9{q) +p,p + a9{q) + p) 



(1) 



where 

d{q) 



1 : g G [0, i) 
— 1 : otherwise 
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with two parameters G M^. The Jacobian matrix of the map reads 



J ^ d<P{q,p) ^1^1- 2a{5{q) + 5{q - \)} 1 
d{q,p) y _2a{5{q) + 5{q-\)} 1 

so that det J = 1 and tr J = 2[1 — a{5{q) + 5{q — \)}]. The map is area-preserving 
and piecewise parabohc. However it is not continuous, the 'discontinuity set' V is 
the following codimension one manifold 

I) = { (0,p) b G [0, 1) } U { I p e [0, 1) } =: v{Q) U v{\) . 

The action of the map can be decomposed into three simple transformations: 

(1) cut of the torus along two vertical lines f (0) and f (|) and translation of the 
two resulting pieces in opposite directions parallel to the cuts: q —>■ q and 

p->p + ae{q); 

(2) rigid translation by (3 along the p-direction: q ^ q and p p + (3; 

(3) parabolic skew translation: q —>■ q + p and p p. 

The triangle map defines a time-discrete dynamical system 0* for t G Z by a recur- 
sion relation 

xt = 0(xt_i) = 0*(xo) , xt=(gt,pt). (2) 
where the t-th iteration of the map can be explicitly written as 

qt{qo,Po)=qo+Pot + - t(t + 1) + a ^ S'fc , 

^ fc=i 

Ptiqo,Po)=Po + /3t + aSt, (3) 
where 

St := d{qk) satisfies 5*4+1 = 5* ± 1 . (4) 

k=0 



Note that 5* G Z and 15**1 < t. 

Using ([3]) and the fact that the function 9{q) is locally constant, it is easy to see 
that the dynamics transforms a horizontal segment to a finite number of horizontal 
segments. On the other hand, for each t, the image of a given vertical line v{qQ) = 
{{qo,p) iP G (0,1]} under the map 0* is a family of segments having a slope 1/t 
and both the length and the position are determined by the pseudo-random series 

i:\-Ak{f3 + ae{qu)). 



3 



Away from the singular set V, the map acts locally as a linear stretching: a 
small parallelogram P, which is bounded by two horizontal sides and two non- 
horizontal sides forming an angle 7 with the horizontal axis, is mapped to another 
parallelogram P'. P' has the same area as P and is again bounded by two horizontal 
sides of unchanged length, whereas the other two sides are rotated clockwise forming 
an angle 7' = 7/(1 + 7) with the horizontal axis, so that they get both stretched by 
a factor sin 7 ■ yT^i-~(r^i-~cot7)^. 

Two nearby points and (g + Ag,p + Ap) on the same side of discontinuity are 
evolving so that 

= ||0*(g + Ag,p + Ap) = ^Ap^ + (Ag + tAp)^ . (5) 



Notice that the distance between points lying on a horizontal segment ( Ap = 0) does 
not grow with time, A(t) = |Ag|, whereas for points on a vertical segment (Ag = 0) 
the distance grows as A(t) = | Ap|-\/l + 1"^- This implies that two arbitrary points 
on the same horizontal segment can be separated only by the cutting mechanism on 
the singular line. This implies that the points in a horizontal line are not separated 
by stretching, but rather by chance that they fall with time on opposite sides of 
discontinuity, which is guaranteed by the stochastic like behaviour which shall be 
discussed in the sequel. 



The triangle map has been further explored and applied in investig ations of fun- 



damen ta l properties of statistic al mechanics in several recent papers (jCasati et al. 



2005bl lal: IPuarte Queiroi . l2008h . In the following text we study, mainly numerically. 



ergodic properties of the tr iangle map using two different ways to symbolically en- 
code the dynamics, see e.g (IChristiansen and Politil . Il996l ). called the polygonal and 
the binary description. We establish some common properties of the two descrip- 
tions as well as some their specific features. In the frame of a given description we 
study a finite state Markov chain corresponding to the triangle map and its spec- 
tral gap. We also analyze certain interesting scaling relations of the spectrum of the 
Koopman operator in the truncated Fourier basis. Additionally we draw parallels 
between the triangle map and its stochastic version called random triangle map. 
The later possesses similar properties as the deterministic version of the triangle 
map, but enables analytical predictions thanks to the possibility of averaging over 
different realizations. 



2 Some further properties of the triangle map 



The ergodic properties of the triangle map system strongly depend on the arithmetic 
properties of the parameters a and j3. In the following we list all the analyzed cases 
and present the results. 
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2.1 The case a = 



If a = the map reduces to the skew translation 4>{q,p) = {q + P + P,P + P)- As 
is well known, for (3 irrati onal the map is (uniquely) er godic but not even weakly 
mixing (see reference e.g. ( Katok and Hasselblattl . 1995 ) and below) 



2.2 The case a and (3 rational 



If a a nd j3 are rational numbers then the dynamics is pseudo-integrable ( iRichens and Berryl . 



198ll ). For example, if a = n/m, and /3 = r/m, for some integers n,m,r such 



that gcd(?T,, m,r) = 1, then the phase space is foliated into invariant curves 
{{q,p + k/m)\q G [0, 1), A; G {0, 1, . . . , m — 1}}. It is easy to see that similar folia- 
tions exist also for arbitrary rational values of the parameters and the dynamics 
restricted on each set of invariant curves is just an interval exchange transformation 
(lET). See Lemma H] below for a more precise statement. 



2.3 The generic case 



The case in which both parameters are non zero and irrational, i.e. 7^ a, /5 G 
with /3 ^ Z + aZ, will be referred to as the generic case. We have the following 

Lemma 1 In the generic case (p has no periodic orbits. 

Proof. According to ([3]), a necessary condition for cf) to have a periodic point {qo,po) 
of period t > 1 is that pt + a St = (mod 1). But this is clearly impossible under 
the assumption of the lemma. □ 

We are now tempted to argue that in the generic case every orbit {qt,Pt)t>o is uni- 
formly distributed (u.d.) mod 1 in T^. This means that for any lattice point (r, s) G 
Z^ \ {0, 0} the sequence {rqt + spt)t>o is u.d. mod 1 in T. Note that if this is the case 
then, arguing for example using the Weyl criterion (see ( Kuipers and Niederreiterl . 



I974J ). Chap. 1.2), the system would be strictly ergodic, i.e. minimal and uniquely 
ergodic. To be more specific, consider a fixed observable f{x, y) of the form f{x, y) = 
^r,s{x,y), where er^s{q,p) '■= * *^ p\ Assuming the uniform distribution of all 
orbits, together with the Weyl criterion, one gets for any 0-invariant measure u: 

1 

= -J2^rAlt:Pt) ^ when n ^ 00 , 




namely i^(er,s) = for all (r, s) G Z^ \ {0, 0}. For an arbitrary / G L^(T^) we write 
/ = J2r,s'^r,sGr,s SO that !/(/) = Cq^q = fj,{f). Siuce a complete proof of the uniform 
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distribution of all orbits is still missing (see however below for a semi-rigorous 
argument) we can only state the following: 

Conjecture 2 In the generic case the system (T^,0, /x) is strictly ergodic. 



We now sketch an argument, half-rigorous, half-heuristic, to support the above 
conjecture. 

The 0-orbit {qt,Pt)t>o is u.d. mod 1 in if for any lattice point (r, s) G \ {0, 0} 
the sequence rqt + spt is u.d. mod 1 i n T. From classical differe n ce th eorems of the 
theory of uniform distribution (see (IKuipers and Niederreitei] . Il974h . Chap. 1.3) 
we know that rqt + spt is u.d. mod 1 iff for every positi ve integer h the sequence 
: = r{qt+h — qt) + s{pt+h — Pt) has this property (see ( Kuipers and Niederreitej . 



u 



1974), p.26). Note that 



{h) 
u\ = r 



h-\ 



+ S 



hpQ + h^{t + 1) + a/i + /i(/i - 1)^ + a ^(/i - k)e[(ik^t) 

^ k=i 

h-1 

hp + aY^ 9{qk+t) 

k=0 



rhpt+i + a rSh + sSh o 0* /? 



r^— — - + sh 



where 



h-i 



Sh:=Y.ih-k)e{qk). 



k=l 



Whence, the time dependence of ul'^^ is that of rhpt^i plus a term which is bounded 
above by (|a| -|- |/3|)(|r|/i(/i — l)/2 + \s\h). Therefore, since pt = po + j3t + aSt we 
have 



lim 



u 



(h) 



r h I 3 + a lim — 

^ t^oo t 



Now, by the ergodic theorem St/t converges a.e. to an integrable function f*{q,p) 
which is 0- invariant and Jj2 f*{q,p) dqdp = and thus Ut^ /t converges a.e. to 
7 = rh{[3 + f*a). Up to now the argument is rigorous. To complete it heuristically, 
one could use Lemma [T] and the invertibility of to argue that cannot have 
eventually periodic orbits and hence that 7 is irrational. Then, stepping backward 
deduce from this that (Mf'*'')t>o is u.d. mod 1 in T and thereby {qt,Pt)t>o is u.d. mod 
1 in T2. 

We have numerically investigated this behaviour and the results are reported in 
Figured], where, in order to simplify the representation, we just recorded the extreme 
points: 



Smax{t) = max 15'^} 

TG[0,i] 



Smm{t) = min 15'^} . 



(6) 
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0.5t + 0.5 



6 e 



10 



12 



14 



Fig. 1. The extreme values of the sum St in the generic case a = 1/e and /3 = (\/5 — l)/2 
at go = 0.1 and po = 0.2. 

It turns out that in the generic case for all tested initial points the sum St yields 
the same behaviour as for the standard symmetric random walks: 



Sr. 



1 



(7) 



The above result is an indication of the stochastic like behaviour of this non hy- 
perbolic dynamical system. We will soon discuss the case (3 = 0, where a different 
asymptotic behaviour of the sums St will emerge. 

We now turn to mixing properties, providing a strong numerical evidence for the 
following 

Conjecture 3 In the generic case the system (T^,0, yu) is mixing. 

In the next experiments we numerically compute both the auto-correlation of a 
given observable f{q,p) as a time average (for a fixed Xq = {qo,Po) G T^) 

1 N-t 

am] = hm — - Y: /(/(xo))/(/+*(a:o)) , (8) 



and the auto-correlation as a phase space average, i.e. 

a(t)[/] = /i(/°0*-/), (9) 



where we consider smooth observable / with zero mean //(/) = /.j.2 f{q, p) dq dp = 0. 
Whenever the system is ergodic the expressions ([8]) and ([9]) are clearly equivalent, 
but since for some values of the parameters ergodicity is still questionable, we sepa- 
rately calculate and compare both quantities (see in particular the next subsection). 
In agreement with both of the above conjectures in the generic case the numerical 
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Fig. 2. The auto-correlations Ct{t) (a) and Cs(t) (b) using observable 
/(x) = sin(27rg) + sin(27rp) at a = 1/e and (3 = (\/5 — l)/2. In (a) we took 
N = 2^^ and averaged the result over m = 4 • 10^ initial points and in (b) we use a 
Simpson integration scheme with Nx = 500, Ny = 5 • 10^. 



resul ts in Figure [2] show the same results as already reported in (ICasati and Prosenl . 
20001 ). i.e. an algebraic decay of correlations with the law 



Ca{t) = 0{t-r 



a = t, s 



(10) 



This indicates that (T^, 0, /i) is strongly mixing with polynomial rate of decay of cor- 
relations. In numeri c al irn plementation of Cg we use a Simpson integration scheme 
( IStoer and Buhrschl . [20oi ) with additional extrapolation between results obtained 
on grids x Ny and 2Nx x 2Ny. Because the observed scaling laws of the map, 
discussed in the next section, are directly related to the geometrical properties of 
the map, it is preferable to consider A^^^ > A^^. Actually, in the following sections 
we will see that the optimal choice is Ny/N^ cx t and oc t, but this is usually 
not possible, because of the large memory consumption. The possible mechanism 
of mixing should be s imilar to that occurring for general polygonal billiards with 
irrational angles (see (jCasati and Prosenl . 1 19991 )). and the model we are facing here 
clearly shows that not only mixing can appear even when hyperbolicity is absent, 
but furthermore it could emerge just by composing two simple piecewise parabolic 
linear maps. 



2.4 The case (3 = 



If P = 0, the map becomes (f>{q,p) = {q + p + a9{q),p + a9{q)), and for rational 
a = n/m we have the following result. 
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Lemma 4 Let (3 = and a = n/m with n, m G N, gcd(n, m) = 1. For po G T set 

m—l 

A{a,po) := IJ h{{po + ka) mod 1) 

k=0 

with hipo) = {{q,Po),(l ^ T} being a 'horizontal foliation' of the phase space. Let 
moreover Q{a,po) = {Qi, ■ ■ ■ ,Qr) be the partition of A{a,po) into disjoint semi- 
intervals whose endpoints form the set A{a,po) fl (j)~^V. The map restricted to 
A{a,po) acts as an interval exchange transformation on the elements of the partition 
Q{a,pQ). Furthermore, ifpo G Q then there exists an integer k = k{a,po) such that 
0^ = id, otherwise acts as an aperiodic lET. 

Proof. The proof is based on the following observation: setting 

O(a,po) = ((Po + Oit) mod l)^^^; > (H) 



we have |0(a,po)| = for all po G [0, 1), so that the set O(a,po) x [0, 1) is made 
of the m 'levels' /i((po + koi) mod 1), A; = 0, . . . , m — 1, i.e. coincides with A{a,pQ). 
Moreover, the set 

I{a,q,po) = {{po + a St{q,p)) mod l)^^^* (12) 



satisfies I{a,q,po) C 0{a,po). The assertion can now be readily checked. □ 
The case a irrational and (3 = represents a situation quite similar to the one 



observed in the free motion of a ball inside a triangular billiard (jCasati and Prosen 



19991 ). and furthermore the most delicate and less explored situation, from both 
the mathematical and the numerical point of view. We start with a simple remark 
concerning periodic orbits and then we turn to ergodicity and decay of correla- 
tions, where new numerical experiments presented here yield evidence of somewhat 
unexpected an behaviour. 

We point out that for a irrational (and /? = 0) we observe the existence of periodic 
orbits. More precisely, according to ([3]), a point {qo,Po) with either po = or po E 
M\Q is periodic of even period t provided the conditions pot + a Sl^o Sk = and 
St = are both satisfied. The set of periodic orbits of length t is characterised by 
the invariant set dt{q,p) = 0, where {qt,Pt) = 4>\qyP) and 



dt{q,p) = J{mm{\q - qt\,l - \q - qt\}y + {mm{\p - pt\ , 1 - \p - Pt\}y ■ (13) 



Formula ([3]) and the fact that the function 6{q) is locally constant immediately 
implies that periodic orbits of given length come in families forming horizontal 
segments. This structure and the statistics of segment lengths does not depends on 
a, but the individual lentgh of the segments does. This is shown in Figure [3l where 
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Fig. 3. The periodic orbits of triangle map at a = \f2 — 1, e ^ (a,b) and /3 = for period 
lengths t indicated in the figure. 

we draw (primitive) periodic orbits for a = \/2 — 1, e~^, up to period t = 10. The 
intervals of periodic orbits are found first by randomly sampling the torus and then 
by searching for nearby zeros of dt{q,p). The number of intervals of periodic orbits is 
growing and their length is decreasing with the increasing period t. Finding compact 
intervals of primitive periodic orbits can be quite time demanding. Therefore we 
instead search for a related object i.e. invariant intervals on 0*, which represent 
the set for all periodic points of length t' such that t'\t. From numerical results, 
shown in Figure H^, we find that the number of invariant intervals A^gx increases 
with increasing t as 



iV, 



fix 



(14) 



In addition to that we check the number of different values of momentum p in 
invariant intervals Pfix- The results are presented in Figure |Dd and they fit well to 
a power law asymptotics 

Pfix X , where C = 1-28 ± 0.02 . (15) 



The number Pgx increases asymptotically slower than A^bx, which indicates that 
invariant lines at a given momentum p are with increasing t cut in more and more 
pieces. Their number is on average given by the ratio N^/ ~ t^~'^ . In addition 
we check the statistics of lengths / of the invariant horizontal intervals at some t. 
We find that lengths are decreasing with increasing t proportionally to t"'', with 
the same C as before. This behaviour is observed using the cumulative distribution 
of rescaled lengths / which does not change with increasing t , see Figure [51 The 
observed scaling reflects the natural simmetry between the diffusion process in the 
momentum and the cutting process in q induced by the discontinuity line. 
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Fig. 4. The number of invariant intervals Nfi^ (a) and the number in invariant intervals 
visited momentum Nf^^ (b) as the function of time t in non-generic case for different values 
of a. 
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Fig. 5. The cumulative distribution of rescaled lengths L = t^'^^/ of invariant intervals 
calculated for a = 1/e and at various times as denoted in the figure. 

Regarding the rate of growth of the sums St-, we find interesting behaviours closely 
related to similar phe nomena tha t one can find already in the case of irrational 
rotations of the circle ( Isolal . 2006 ). In particular, we find evidence that for a being 
a quadratic irrational, the extreme values of the sum Si (Smax Slid — ^min) increase 
logarithmically, see Figure O We note that to obtain these res ults we had to rely 
to hi gh (more than double) precision fioating point arithmetic (IDemmel and Hida . 
20031 ). Moreover, similarly to the case of rotations discussed in ( llsolal . l2006l ). we 
expect a significant dependence of the growth behaviour of St on the arithmetic 
type of a. 



On th e basis of preliminary numerical experiments performed in (ICasati and Prosed . 
2000l ). it was conjectured that for a irrational and /3 = the system (T^,0, /i) was 
just ergodic, and perhaps only weakly-mixing with no decay of correlations. Fur- 
ther and more detailed numerical experiments reported here indicate that certain 
polynomial decay of correlations can instead be present. This behaviour strongly 
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Fig. 6. The points of extreme values S'ma? 
a = \/2 — 1, (-v/5 — l)/2 (a,b) at go = 0.1 and po 
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Fig. 7. The auto-correlation Ct(t) using observable f{q,p) = sin(27rg) +sin(27rp) at /3 = 
and a = \pl — 1, (a,b), where we take = 2^^ and average the result over m = 10^ 
realisations. 

depends on the values of a, and probably also on the specific choice of the observ- 
able but nevertheless it appears as a solid evidence of mixing. This is yet another 
indication of the multifaced behaviour of the maplH More precisely, in most cases of 
irrational a we find numerically non-exponential mixing with the auto-correlation 
approximately described by 



C{t) = Oit-"^) 



where 00 = 0.2 ± 0.02 . 



(16) 



In Figure [7] we show the auto-correlation Ct{t) for two representative examples 
a = — 1, using observable f{q,p) = sin(27rg) + sin(27rp). We observe a very 
smooth decay. Basically the same results are obtained also for a = {^/E — l)/2 and 
by using other observables e.g. f{q,p) = sin(27rg) + sin(27rj9) + sin(27rg) sin(27rp), 
or f{q,p) = 0{q). The auto-correlation Cs(t) is not shown as it is numerically in- 
distinguishable from Ctit). We also found an exception to the given approximation 
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Fig. 8. The auto-correlations Ct{t) (a) and Cs{t) (b) of the observable 
/(x) = sin(27rg) + sin(27rp) at /3 = and a = \/3 - 1. In (a) we take N = 2^° 
and average the results over m = 4 • 10^ realisations, and in (b) we use the Simpson 
integration scheme with = 200, Ny = 6 • 10^ 

e.g. for a = ^/S — 1, where we find for certain observables, e.g. sin(27rg) + sin(27rp) 
and sin(27rg) + sin(27rp) + sin(27rg) sin(27rp), that the correlations could decay a bit 
faster, approximately as 



c{t) = o(r) 



where = 0.4 ± 0.05 , 



(17) 



whereas for other observables (e.g step function) the law of decay is (fT6|) . The 
anomalous situation in presented in Figure [H where we show results for both defini- 
tions of auto-correlations Ct{t) and Cs{t). It is not surprising that in this situation 
the connection between a, the choice of the observable and the power of the cor- 
relation decay u can be quite complicated. Notice that due to a slow diffusion in 
momentum space, the calculated correlation decay in the presented time windows 
is actually taking place on a small subset of the phase space i.e. possibly on only 
few horizontal intervals per initial condition. 



3 The Koopman operator 



Some further particular cases can be discussed in terms of the spectral properties 
of the unitary operator U on H = L^(T^, dfi), where dfi = dqdp, defined by 

{Uf){q,p) = fi<Piq,p)). (18) 



Consider the orthogonal decomposition H = Hi (B H2, where H2 is the subspace of 
functions depending only on the second coordinate, i.e. the subspace generated by 
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{eo.Jsez where e^,,(g,p) = e^^^^'-'^+^p). By the relation 

U' er,s = er,s+rt C^^' Ello ( ° , (19) 

which for t = 1 writes 

Uer,s =e,,,+,e2'^'(^+'-)("^(^)+'^). (20) 

On the other hand, a simple calculation yields 

g2.i.(.%)+/3)^g2.i./3 / (27r£a) + -sin(27r£a) V ^] (21) 

V ^ fc odd ^ / 

and therefore 



Uer,s =e'"'(^+'^)^ (^cos [27r(s + r)a] e 



r,s+r 



+ -sm[27r{s + r)a] ^ ) (22) 

^ fc odd ^ ^ 

From this expression we immediately deduce the following 

Lemma 5 If a = n/m and /3 G M then U has eigenfunction eo,m; with eigenvalue 

Proof. Writing fl22l) for a = ^, r = and s = m we get U eo,m = e^'^'™^ eo,m- ^ 

Corollary 6 Let /i be the two dimensional Lebesgue measure. If a = n/m and 
l3 E M. then (¥^,0, yu) is not weakly-mixing. If furthermore /3 G Z + aZ then the 
system is non-ergodic. 

Now, given a function / G L^(T^) let us expand it as 

f = ^Cr,ser,s with ^|Cr,sP< 00 (23) 



so that by fl22l) we have 



f/ / = E Cr,. e2-i(^+'^)/5 COS [27r(s + r)a] e,,,+. 

+ ^sin[27r(s + r)a] ^ ) . (24) 

A; odd 
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Suppose now that there is A 7^ with |A| = 1 such that Uf = A /. By ( !23|) and ([2 
the eigenvalue equation becomes 



Ac^,, = e^^'"^ (cos(27rsa) c^,,_^ + - sin(27rsa) ' ) ■ (25) 

\ ' vr fc odd ^ 



For a = equation fl2^ becomes 

Thus, when restricted to the (invariant) space H2-, U multiphes the Fourier coeffi- 
cients by e^'^'^'^ and therefore has a discrete spectrum. Instead, on Hi (where r 7^ 0) 
U generates an infinite family of trajectories with infinite cardinality and therefore 
has a Lebesgue spectrum. In summary, for a = and (3 irrational (F, T^,/i) is er- 
godic but not weak mixing (althoug h for functions belonging to if i one may have 
arbitrarily fast decay of correlations ( Courbage and Hamdan . 19971 )). 



For a 7^ the space H2 is not invariant. An eigenfunction / G would give 

Aco,s = e^""'"^ cos(27rsa) co,s 

which is impossible for a irrational, since it would be |A| < 1 (unless / = Co,o 7^ 

and A = 1). Moreover, for s = we get A Cy Q Cj" — f SO that Icr^ol = |cr,-r|- Siuce 

Sr,s |cr,sP < 00 we have Cr,o = Cr-r = for all r 7^ 0. Using Cr,o = and fl25|) with 
r = s we see that Cr^r = for all r 7^ as well. Therefore, if co,o = (so that A 7^ 1) 
then / must be of the form 

|r|^|.|^0 



where the sum has at least one term with r 7^ 0. 



3. 1 Fidelity of the truncated dynamics 



Having no further analytical results concerning the spectral properties of the op- 
erator [/, we now present some interesting scaling laws observed in the numerical 
investigation of a suitable truncated Koopman operator Utr on the finite dimensional 
space i^tr spanned by the Fourier basis {er,s}\r\<N,\s\<M- Even if at the moment we 
do not have a clever theoretical understanding of these results, their evidence is so 
robust that we believe it is worth mentioning. 

The operator f/tr is no longer unitary. The Fourier basis is chosen so that f/tr? as 
previously f/, maps real functions into real functions satisfying the identity 

KU,,K = Ut, , (28) 
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where K is the operator of the complex conjugation. Thereby complex eigenvalues 
of Utr come in complex conjugated pairs. The parameters and M determine the 
dimension of H^^ and we express them by using new parameters D and k as 

N = {{2D + - 1) /2 and M = {{2D + 1)^-^ - l)/2 , (29) 

so that {2D + 1)2 = {2N + 1)(2M + 1) and k = {2N + 1)/(2M + 1). We see that 
D represents the geometric mean and k, the ratio between the numbers of Fourier 
modes chosen along q and p axis. We quantify the fidelity of the truncated dynamics 
generated by the Utr by 

Fn{t) = \Kgh , F^{t) = \\g o 0* - Ulgh , (30) 

applied to an observable g (and its truncation in iftr), with Jj2dqdpg{q,p) = 0. 
The first quantity, -Fu(^); measures the deviation from unitarity, whereas the second, 
-Fd(t), captures the distance from the dynamics generated by the triangle map. Here 
we use an observable g{q,p), which is constructed by projecting a two dimensional 
Gaussian function, centered in {qo,Po) and of variance s~^, onto the torus and 
subtracting its mean. Thereby we obtain 

9{(1,P) = 9s{q-qo)9s{p-Po) -'^ , gs{x) = J-J2exp(^-s{x-nf) . (31) 

' ^ nez 

The L2 norm of the Gaussian packet is ||5'(q',p)||2 = 6'jac(27r2/s)2, where 6'jac(a;) = 
J2nez exp(— X ri^) is connected to the Jacobi theta functions. Nevertheless, in the fol- 
lowing, presented results are similar for all sufficiently smooth functions i.e. those 
with the Fourier expansion decaying sufficiently fast. In Figure [9] we show the time 
evolution of -Fu(^) for one Gaussian packet. The analysis of results reveals a remark- 
able scaling property of F^^{t) w.r.t. the parameters D and k valid in a certain time 
window <^ t < 0{\/15), and reads 

1 

log(Fu(t)) = -a^t + b for D > 1 . (32) 

where a and b are constants (independent of t,D,K). Surprisingly, the constant 
a is not strongly dependent on observables and its value is around | for various 
observables that have been tested. The obtained scaling indicates that the property 
of unitarity decays exponentially approximately with rescaled time t/dim{Ht,.)i « 
t/ D^/"^ ~ and the conservation of the norm is improving with decreasing k. A similar 
scaling law is observed in the quantity -Fd(^); which at t = starts from zero and 
then with increasing time saturates to a plateau given by H^? o 0*||2 = ||(?||2, see 
Figure [10 



16 




Fig. 9. The evolution of F^^(t) for a Gaussian packet g {s = 10, go = 0.3, = 0.6) at 
a = 1/e, P = (\/5 — l)/2, where in (a) it is calcuated at k = 1 and for different D (a), 
and in (b) at fixed D = 500 and for various k (b) as denoted in the figure. The fitted 
constants are a = 0.549624 ± 0.0008 and b = 0.578 ± 0.007. 
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Fig. 10. The L2 norm of the difference between a Gaussian packet g 
(s = 10, go = 0.3, po = 0.6) propagated by the truncated Koopman operator and 
by exact triangle map calculated for a = 1/e, /? = (\/5 — l)/2 with k = 1 and for various 
D. 

3.2 The spectrum of the truncated Koopman operator 




The spectrum {A : det(f/tr — Aid) = 0} of the truncated Koopman operator on 
this finite functional space of dimension {2D + 1)^ is inside the unit-circle on the 
complex plane. The side effect of the truncation is that Utr is not invertible, having 
the rank smaller than the dimension, see Figure [TTb . Notice that only non-zero 
eigenvalues are far enough from the unit circle participate in decay of F^^{t) from 
the start. We studied the distribution of modulus of eigenvalues denoted by px{x) 
and have found that the probability scales like 0{D~2)^ see Figure [TTb. On the 
other hand the distribution of eigenvalue modules has singularity at zero of the form 
P\{x) = 0{x~^). It is worthwhile to mention that eigenvalues have higher density 
near the unit circle. The peak of eigenvalues is with increasing D slowly moving to 
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Fig. 11. The rank of the truncated Koopman operator as a function of D (a) and the 
distribution of the modulus of eigenvalues px{t) for different D (b) calculated at a = 1/e 
and (5 = (\/5 - l)/2. 
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Fig. 12. The largest eigenvalues of the truncated Koopman operator for different D (a) 
and the number of eigenvalues larger then some p as a function of D (b) both calculated 
at Q = 1/e and (3 = (^/5 - l)/2. 



the unit-circle and thereby we observe the squeezing of the spectral gap defined as 
A = 1 — max{|A| < 1}. The gap in our case is a non-smooth function of D. However, 
due to the condensation, the gap does not have an important dynamical meaning. 
The number of eigenvalues by magnitude larger than some value < x <C 1 is 
increasing as 0{D'^), where fitted value of the exponent is z/ = 1.634 ± 0.002, see 
Figure fT2k . This scaling can be also applied to eigenvalues near the unit-circle so 
that behaviour of eigenvalue modules near the circle align for different D as we 
can see in Figure [T2b. Although the connection between the fidelity i^u(^) and the 
spectrum of the truncated Koopman operator is clear, we find that the behaviour 
of eigenvalue modulus with D is insufficient to explain the scaling properties of the 
fidelity. 
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4 Dynamics of polygons 



In this section we consider a heuristic but quite fruitful geometric approach, namely 
we study the dynamics of the map acting on particular simple subsets of phase 
space ^ C T^, the polygons. It is easy to verify that all the ingredients of the 
dynamics preserve polygons, namely maps a polygon to a polygon, and hence 
its iterates 0" map a polygon to a finite set of several polygons. Furthermore it 
is possible to make an efficient and accurate numerical algorithm for performing 
dynamics of polygons by representing them - or their boundaries - in terms of sets 
of lines. 

In the following subsection we describe and apply such explicit construction and 
the related numerical procedure in order to study the statistical properties of trans- 
formed polygons as it emerges from iterating the dynamics for sufficiently long 
times. 

4-1 Statistical properties of polygons 

We now consider an arbitrary initial polygon ^ C and evolve it with time 
At := 0*(^), t G Z. After some time t, the set At is composed of many, say N(t), 
pieces. The following main questions are addressed in our numerical experiment: 

• Do there exist well defined statistical distributions of geometric properties of 
these little polygonal pieces, such as area, diameter in x-direction, diameter in 
^/-direction etc., after long time t, and do these distributions possess properly 
defined averages and variances? 

• What are the scalings of geometric properties of little pieces as t grows? Are 
asymptotic, t —>■ oo, statistical distributions universal, i.e. independent of the 
geometry and volume of the initial set A7 

• Are these pieces distributed uniformly in phase space T^, i.e. that the number 
of pieces inside an arbitrary ('window') set B G T"^ with area fJ^iB) is, after long 
time t, equal to ^{B)N{t)l 

A computer program has been written which computes an exact evolution of an 
arbitrary initial polygon A in terms of a set of polygons At- Affirmative answers to 
the above questions (for generic, irrational values of parameters a, (3) are strongly 
suggested by the numerical results which are summarised below. It can be easily 
shown that any connected piece, say P, of At has the following properties (which 
are preserved under time evolution) : 

^ Provided such conditions are imposed on the initial polygon = A. Otherwise, these 
properties hold for ever increasing fraction of its little pieces, namely for those whose 
boundary no longer contains the image of any part of the boundary of the initial set A. 



19 



(1) There exist left-lower and right-upper corner points, Ci,C2 G V, with coor- 
dinates {qi,Pi) and (52,^2), respectively, such that for arbitrary {q,p) G V: 

qi < q < q2, Pi <p < P2- 

(2) The boundary of the polygon V can be written in terms of two piece-wise 
linear curves (monotonically increasing in q — p plane) going from Ci to C2, 
namely lower B\ and upper boundary. Going from Ci to C2, the slopes of the 
straight segments of in g — p plane can be written in terms of an increasing 
sequence of positive integers 1 < rii < n2 . . ., as l/rii > l/n2 > ■ ■ ., while the 
slopes of the straight lines of Bi are given by a decreasing sequence of positive 
integers rrii > m2 . . . > 1, as 1/mi < l/m2 < . . . 

In Figureliniwe show the image of the initial polygon (first plot) in terms of five snap- 
shots (additional five plots) at generic irrational values, a = e~^,P = (a/5 — l)/2, 
where integer times t shown are roughly by a constant factor apart. Each snapshot 
is plotted in a renormalised window, around certain reference point {qc,Pc) = (1, 1) 
(in the center of the phase space, which is here scaled to [0, 2) x [0, 2)), of g— width 
6g/t and p— width 6p/t'^ in order to suggest the scaling. 

Observe the statistical similarity of consecutive snapshots! 

By repeating the experiment for a smaller initial set the statistical self-similarity 
persists just the density of little pieces reduces proportionally to the area of the 
initial set. Figure [H] presents computation for the non-generic ('weakly-erg odic') 



case of (3 = (such map has been studied already in (iKaplan and Hellerl . ll998l )). We 
see large temporal and spatial fluctuations in the density of little pieces which may 
complicate the mechanism of dynamical mixing in this situation. Clearly, general 
heuristic arguments for the mixing mechanism would not hold for this case. In the 
other non-generic case a G Q and (3 G M\Q we get similar nonuniform structure 
with disappearance of some stripes in initial transient time. 

We count the number of the little pieces N(t) obtained by evolving initial set up to 
time t. In Fig{T5]we show N{t) for the two generic cases with different initial sets and 
for a weakly ergodic case. Always, the growth of N{t) is found to be cubic in time 
after some transient period, with the proportionality constant factor proportional 
to the area of the initial set 

N{t) ^Ct^ , C(x ^i{A) . (33) 



We analyse the statistical distributions of the following geometric measures of little 
pieces: (i) the area a := fi(V), (ii) the g-diameter Sg := q2 — qi, (iii) the p-diameter 
'■= P2 — Pi, and (iv) the average slope of the piece s := Sp/6g. 

(a) =30.0/t^ aa = 1.12(a) , (34) 
{6g) = 7.2/t , as, = 0.59 (5,) , (35) 
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t=284 t=455 



Fig. 13. Renormalised (zoomed) snapshots around some point in phase space 

{qc = ^,Pc = 1) at different times on a rcscalcd total phase space [0,2) x [0,2), for 
generic case a = e~^,P = {y/5 — l)/2. Vertical labelling to the left of each plot indicates 
the coordinates {{qi,q2},{pi,P2}} of the window shown. 
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Fig. 14. Renormalised snapshots in non-generic case a = e , /3 = 0. See also caption of 

{5p) = lAA/e , as^ = 0.52 {5^) , (36) 
(s) = 2.2/t, a, = 0.37 (s) . (37) 

This quantitative results and self-similarity of pieces on the phase space suggest a 
kind of statistical scaling invariance under the transformation q' = qt and p' = pt^. 
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Fig. 15. The number of pieces A^(t) as a function of iteration-time t for three different 
cases of initial set sizes and their estimated asymptiotics Ct^. 



Therefore, let us introduce t-invariant quantities 



Xq 



t6, 



Xp 



Xs 



ts . 



(38) 



Next we investigate the full distributions of scaled quantities Xa,q,p,s- 
by calculating two different distributions: 



This is done 



(i) Cumulative number distributions, equation for the quantity Q which is one of 

Xa? Xxj Xj/; Xs 



Nq{x) 



#{g(p) < x} 

N{t) 



(39) 



(ii) Cumulative probability distributions, 

where H{x) = {0 : a; < 0; 1 : x > 0} is the step function. 



(40) 



We note that the qualitative features of the number distributions Nq and probabil- 
ity (measure) distributions Wq are essentially the same, so we shall in the following 
only report numerical results on the former. It has been checked with great numer- 
ical accuracy that distributions Nq{x) and Wq{x) are asymptotically (for large t) 
independent: (i) of time t, (ii) of the size and shape of the initial set A and (iii) of the 
particular values of generic irrationals a, (3. In Figure Uni we show numerical results 
on these distributions which suggest several small and large argument asymptotics 
which are indicated within figures' labels. It should be noted that the numerical 
data in Figure [16] can be described even globally quite well by an exponential fit 



Nx,y{x) = 1 - exp(-7a;,j/X) with some exponents 7a;, 
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Fig. 16. The cumulative number distribution of the polygon areas Na{x,t) (a), polygon 
horizontal dimension Nx{x,t) (b), and polygon vertical dimension Ny{x,t) (c), for dif- 
ferent times t in the generic case at a = 1/e and f] = - l)/2. The fitted curve is 
f{x) = -1.6644 X + 0.475028. 



We conjecture and show later in subsection 16.21 for a simplified special case that 
these distributions (and the scaling) can be derived from the fixed point condition 
for certain dynamic renormalisation group equations. The later is found for the 
random triangle model or more precisely for its one-dimensional simplification. 
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4-2 A Markov approximation 



Let us define a finite Markov approximation of the dynamics by the following con- 
struction of a transition matrix. We start by the set of two vertical discontinuity 
lines T>. They can be viewed as boundaries of two elementary polygons - vertical 
stripes - which divide the entire torus T in two pieces and have no corners. This 
is just the initial partition that we will be using later to construct a suitable bi- 
nary coding of the dynamics. Here we discuss only the generic case at parameters 
« = 1/e, /3= (V5- l)/2. 

Let us now study subsequent images of this set and their union up to integer time 
t: 

t-i 

= (J (41) 

n=0 



A set of lines V^*^ can be considered as a boundary set for a finite set of N(t) 
polygons, let us denote it by Al^*^. It is obvious that : M^^^ M^^'^^^ . For 
any t, a set A4^^^ defines a convenient partition of phase space to N{t) elements 
on which we define a probabilistic process in terms of a Markov transition matrix 
M(*) of dimension N{t). Note that C if t' < t. It is obvious that each 
element A & A4^*^ is mapped again onto a single element B e A^^*'', provided this 
element is bounded only by the lines from the set I>^* ^ for t' < t. For each such 
A the column of transition matrix reads M^,-*^ = 5^' ,5. Hence M^*) there acts as 
a simple permutation. Note that since the number of elements grows as N{t) oc t^, 
this happens for a majority of cases, namely only for N'(t) oc N{t + 1) — N{t) oc 
elements A E A4^*^ it may happen that their boundary also includes lines from the 
last image (f>^*\V) and hence they are cut and split into several elements of M.^*\ 
Since the boundary set T>^*'> is a subset of , it is obvious that each such polygon 
A can be cut only once, and hence mapped to at most two polygons, say Bi and 
B2. The probabilities are given simply in terms of area measures of image polygons, 
hence such columns A of the transition matrix read 



-U'-^JdxldyS{y-<l>{x)) 



A A' 

- ^ {Sa',bM^i) + Sa',bM^2)) . (42) 
To each split of a polygon we associate the relative sphtting strength as 

ri(A) := mm{i,(B,),ii(B2)}/li{A) < \ . (43) 
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Fig. 17. The distribution of free propagation length of partitions (a) and the cumulative 
distribution of splitting strengths (b) in the Markov matrices calculated in the polygonal 
description at a = 1/e, /? = (\/5 — l)/2 and for different times t. 

We study the cumulative distribution of splitting strengths, denoted by Prf{x). In 
Figure [T7k we plot the Pri{x) for different t. We notice that small splitting strengths 
are dominant and the probability density distribution has a square-root singularity 
at 7] = 0. The distribution converges with increasing t to its limiting form, which is 
well fitted by the phenomenological formula 

Nr,{x)=y + ay{l-y){y + b)+e, y = , (44) 



2 

0.2 0.4 0.6 0.8 1 

(b) 



with absolute error |e| < 5 ■ 10 ^ and constants a ^ 0.1 ± 0.05 and b ^ 2± 0.3. The 
first two central moments of the limiting distribution are 

7/ = 0.2 ± 0.002 , = 0.155 ±0.004. (45) 



In the following we are interested in dynamical and statistical properties of the 
Markov process for finite times t and how they evolve with time t. The Markov chain 
is found to be ergodic and mixing with the spectrum denoted by {Xi{t) G C}fJ;Q~^ . 
In Figure [18] we show two examples of spectra of finite Markov matrices and the 
distribution of eigenvalue magnitudes |Aj| which shall be later compared to Markov 
spectra with respect to an alternative partition coding. Further on, we study the 
evolution of the spectral gap 

A(t) = 1-1 max{Ai(t) ^ 1}| (46) 

i 

and Kolmogorov- Sinai entropy 

h{t) = - y: m e <^^'iogK!^')- (47) 
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Fig. 18. The eigenvalues Aj of the Markov matrix in the complex plane calculated in 
polygonal description at a = 1/e and P = (\/5 — l)/2 for t = 15 (top, left),t = 30 (top, 
right) and the distribution of eigenvalue amplitudes px{x) in different representations 
(bottom). 

Because the triangle map is non-hyperbolic, we expect that in the limit t ^ oo the 
Markov matrix does not have a spectral gap and entropy is equal to zero. This is 
supported by the numerical results shown in Figure [191 The entropy as a function 
of time decreases monotonically. Its asymptotic dependence fits well to an algebraic 
law 



h 



~ C'entropy'^ 



where r = 0.97 ± 0.03 . 



(4^ 



Something similar can not be said for the gap, because its dependence on time is 
not so simple. Nevertheless numerical results suggest that the gap asymptotically 
decreases inversely proportional with time as 

A~Cgapt"'. (49) 
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Fig. 19. The largest amplitude eigenvalue for different times (a) and the spectral gap A 
together with the Kolmogorov-Sinai entropy h as a function of time (b) in the Markov 
matrix in polygonal description at a = 1/e and /3 = (\/5 — l)/2. The fit to the entropy 
curve is equal to /ifit(t) = 2.529 1-^-^^^. 

5 A suitable binary coding 



The description in terms of polygons leaves a lot of freedom in the choice of the 
initial setup. Nevertheless, the dynamical properties of time asymptotics are shown 
to be only trivially dependent on the initial partition set. In the following we 
present a description of dynamics base d on binary coding of t rajectories. In cer- 



tain chaot ic systems, e.g. baker m ap (lArnold and Avez 



map (see (ICvitanovic et al.l . Il988l ). (jP'Alessandro et al 



199lh ) 



Lozi and Henon 
such description 

is proved to be asymptotically exact. We shall provide evidence that this is the case 
for our non-hyperbolic system as well. 



5. 1 Topological complexity of the triangle map 



Let A = {Aq, Al} be the binary partition of into the half-open rectangles defined 

as 



A, = [^,^) X [0,1), z = 0,l 



(50) 



and we consider the natural coding associated to it. More formally, we let a : 
{0, 1}^ be the coding map defined by 



(51) 
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where Xa{(1,p) = (1 '■ {liV) G : otherwise) is the characteristic function of set 
A. One may wonder if A is a generating partition, in the sense that 

cT(g,p) = a{q',p') {q,p) = {q',p') . 

Since horizontal segments remain horizontal during the evolution, whereas arbitrary 
small vertical segments are linearly stretched while they become asymptotically 
horizontal, one easily checks that 

a{q,p) = a{q',p')^p = p' . (52) 



Therefore, in order that A be generating it suffices that the sequence {qt)t>o defined 
in ([3]) is dense in T for all {q,p) G T^. In the non-generic case, a = and /3 G M\Q, 
this property follows at once from a theorem of Weyl according to which if P{x) is a 
polynomial of degree m > 1 with real coefficients and at least one of the coefficients 
is irrational then the sequence P{n) is uniformly distributed mod 1 and thus dense, 
whereas in the generic case this is directly related to the validity of Conjecture [2j 

We associate to an individual binary code uj = {uq, ■ ■ ■ ,ujt-i) G {0,1}* a set of 
points A^ defined as 

t-i 

A^ = {x : (P^ix) e A,„ < A; < t} = fl 0"' (Ac.J • (53) 

A;=0 



We denote by Lf C {0, 1}* the set of 0-admissible words of length t, i.e. 

Lt = {uje {0, 1}* : A^ ^ 0} . (54) 



The partitions A*^*) based on the binary partition A is then defined as a collection 
of all non-empty A^ 

A^') = {A.LeL. • (55) 



Note that each is a union of disjoint polygons with horizontal parallel sides, 
whereas the other two opposite sides are in general not parallel and each has a 
slope not smaller than 1/t (see also below). In Figure [201 snapshots of the iterated 
partitions A*^*^ for different values of t are shown, in fact we plot the set of polyg- 
onal boundaries of all elements of A*^*\ which is just P'-*^ (see equation (HT]) ). The 
cardinality of the set Lt, that is 

N{t) = if{iu e {0, 1}* : A^ ^ 0} , (56) 
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Fig. 20. In the upper two plots we show polygonal boundaries of the partitions A^*^ for two 
subsequent time steps, namely P'^*-' (blue thick lines) and on top of them (red thin 

lines), for t = 3 (left plot) and t = 4 (right plot). Fundamental square [0, 1)^ is indicated 
by a thin black line, and periodic boundary conditions are imposed. In the lower plot we 
show a finer partition, i.e. polygonal boundaries 'D^^\ 

is called the topological complexity function of the pair (0, A). In positive entropy 
systems the cardinality of the set is expected to grow exponentially as 

N{t) = 0{exp{ht)) , (57) 



where h is the topological entropy. In the non-hyperbolic triangle map we are con- 
sidering here h = and therefore we could expect only an algebraic growth. An 
interesting problem is that of estimating the asymptotic growth of N{t) as t — >■ oo 
which is up to a constant factor equal to N{t) discussed in Section HI The numer- 
ical experiments (see Figure El]) lead to the conjecture on existence of a positive 
constant C independent of the value of irrational parameters {a, (3) such that 

N{t) r^Ct^ , C = 0.66 ± 0.02 . (58) 



A slightly weaker result was recently proved in (IBonannd . l2006l ) and argues that for 
almost all a and /3 one can find two positive constants Q, i = 1, 2 such that Ci t^ < 
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Fig. 21. The number of partition elements N(t) for various values of parameters a,P 
as indicated in labels calculated by counting the binary codes corresponding to points 
uniformly sampled over the torus. 

N(t) < C'^ t^. Moreove r, a similar result is known in rational polygonal billiards 
(ICassaigne et al.l . |2002| ). Compatible with this law, at time t about elements of 
the partition get cut into two new pieces, i.e. for about words {ojq, . . . ,ujt-i) we 
have 



A, 



U)(,,...,U)„--l 



(59) 



with both Aa,o,...,c^„_i,o 7^ and Ai^o,...,a;„_i,i 7^ 0- 

The partition elements A G A^*) are of different sizes (measures) /i(v4) and are 
decomposed of disjoint polygons, see Figure [22b. We calculate the (heuristic) prob- 
ability Pm for a partition element being composed of exactly m G N disjoint poly- 
gons. An example of the latter is presented in Figure [^2b for the generic case. The 
typical number of polygons in a partition element is one or two, where the second 
is slightly more favoured. In addition we see that the probability Pm is decreasing 
approximately stretch-exponentially, logPm = ~\0{\/rn)\. The maximal number of 
polygons in a partition element appears to be growing approximately proportionally 
with t, but this result is not strongly statistically relevant. The presented properties 
about decomposition of partitions for the generic case do not change significantly in 
the non-generic case. The areas of partition elements are decreasing with increasing 
t. We study their statistical properties using the cumulative distribution of areas 

iVa(x, t) = i^{ueLt : x> /i(A^)}/iV(t) . (60) 



and cumulative distribution of rescaled areas Np{x,t) = Na{N(t)x,t). In Figure [231 
we present Np{x, t) in the generic case for different t. We see that Np{x, t) converges 
with increasing t to an unique distribution A^p(a;) with the following behaviour 

Np{x) = 0{x^) for X < 1 , (61) 
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Fig. 22. The elements of the partition of the triangle map at time t = 4 individually 
coloured (a) and the distribution Pm of the number of polygons m in the elements of 
the partition (b) obtained at different times t indicated in the figure, both obtained at 
a = 1/e, p = (Vs - l)/2. 




Fig. 23. The cumulative distribution of rescaled areas of partition elements Np(x,t) in 
triangle map in two different representations calculated at a = 1/e and /? = (\/5 — l)/2 
for different length of trajectories t. The fitted curve is f{x) = -1.88098 x + 0.945298. 

and 

Np{x) = e-l^("^)l for x > 1 . (62) 



Consequently the average area of partition elements in average decreases as 0(1/ N(t)) 
with increasing t, which is similar to the polygonal description. 



5.2 Transition probability matrix 



We are now again going to approximate the dynamics of the triangle map through 
a probabilistic model i.e Markov chain. The long range correlations in this system 
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prevent (statistical) validity of this description on the long time scale. Nevertheless 
it gives an useful insight into the properties of the dynamics and the slow mixing 
decay which we numerically detected in the triangle map. 



Simil ar as above fl42p we associate to the triangle map a Markov process (IMeyn and Tweedid . 



19951 ) build on transition probabilities 

M,,,, = Prob (0(x G A^) G A^O (63) 



between elements of the partition A*^*^ = {A^ : cu E Lt}. One has 
1 

MA 



<L' = I dx / dy6{y - 0(x)) . (64) 



Aa 



The transition probabilities determine the Markov matrix M*^*) of our process. Be- 
cause the cardinality of the partition N{t) grows as 0{t^) with increasing t, the ma- 
trix M*^*) is sparse for large t and therefore can be efficiently stored. The sparseness 
is increasing with increasing t. We recognise that the number of non-zero elements 
in M^*^ is equal to 

#{<L' ^0: u,uj' ELt} = Nit + 1) . (65) 



By considering the Markov matrix of the mth iterate of the map the number of 
non-zero elements is N(t + m). Note that a point x G A^^ with u = {uq, ui, . . . , Ut-i) 
is mapped to 0(a;) G A^^/, where a;' = {uJi,U2, ■ ■ ■ ,ujt-i,€-) and e G {0,1}. Conse- 
quently each column of M*^*^ has at most two nonzero entries denoted by r]^ G (0, |) 
and 1 — T]^. We refer to rj^^ as the splitting strength of the element in A^*) labeled by 
UJ and is defined as 

7]^ := min M^^l, . (66) 



An individual element of the partition A*^*) is mapped by the Markov matrix M^^^ 
usually to another element or its image is split between two elements. We calculate 
for each element in A*^*-* the free propagation length denoted by /(e) , e G A*^*) i.e. how 
long it can be evolved by M*^*) without being split. This time is necessarily finite as 
the opposite would indicate that the map is non-ergodic. The distributions of free 
propagation lengths it is defined as 

Pi{x,t) = > (67) 



and shown in Figure [2lk for different t. Properly scaled free propagation lengths 
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Fig. 24. Distribution of free propagation lengths of the binary partition A*^*-* (a) and 
cumulative distribution of splitting strength Prj{x,t) of the Markov matrix (b) calculated 
at a = 1/e and /? = (\/5 — l)/2 and different times t as indicated in the figure. 

reveal that the distribution Pi{x,t) has an asymptotic form for large t 

P,(a:,t) = ^P.(|) , Pz(x) = e-l«(-')l , (68) 

We see that the topology of transitions in the Markov matrix strongly suppresses 
a long free propagation. The average free length scales with time as (/) ~ const t, 
which can be read from fl68l) . The splitting strength strongly varies between parti- 
tions. This can be seen in Figure [2Dd, where we plot the cumulative distributions 
of splitting strength 



for one set of parameters. We see that the variable is almost uniformly dis- 
tributed and that the distribution of r] converges with increasing t to a form well 
described by model (gl]) with a = 0.195 ± 0.005 and b = 1.42 ± 0.05. It seems that 
the limiting distribution is common for all generic cases with the first two central 
moments reading 

rj = 0.2 ± 0.005 , a^ = 0.15 ±0.01. (70) 



Without detailed analysis we could say that the the distribution of splitting strengths 
in the binary and polygonal description is almost identical. 

The eigenvalue spectrum E(t) = {A G C : det(M*^*^ — Aid) = 0} of the Markov 
matrix is confined to the unit circle. Particularly interesting is the spectral gap 
A between the unit-circle and the second largest eigenvalue. Examples of spectra 
obtained for different times are shown in Figure [251 We notice a persistent eigen- 
value A = 1 corresponding to the invariant distribution. Other eigenvalues inside the 
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Fig. 25. The eigenvalues Aj of the Markov matrix calculated using the binary description 
at a = 1/e, P = (\/5 — l)/2 for times t = 16,30 (top) and the distribution of eigenvalue 
amplitudes in semi-log and log- log presentation (bottom). 

unit-circle move with increasing time t toward the unit circle. This can be clearly 
observed in distribution of eigenvalue amplitudes 

p,{x,t) = -^ Y: 6ix-\X\), (71) 
y^) AeE{t) 



and in particular by the spectral gap A{t) as a function of time t. The general 
form of the distribution of eigenvalues amplitudes do not change with parameters a 
and (3. An example is shown in [251 The concentration of eigenvalues increases near 
to the unit circle and around the origin of the complex plane. Around the origin 
the distribution has an algebraic singularity, approximately px{x,t) ~ const x~^, 
similarly as in the case of the Koopman operator. From the Figure [2Bb we can 
read that the compact "cloud" of largest eigenvalues in average converges to the 
unit circle as 0(t~^). The outermost point defines the spectral gap of the system, 
which is plotted as a function of time in Figure [26b. We see that the gap decreases 
with increasing t algebraically, nevertheless its functional dependence on time t is 
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Fig. 26. The largest eigenvalues (a) and the spectral gap with the entropy of the Markov 
matrix (b) at parameters a = 1/e, /? = (-v/S — l)/2. The fitted curve correspond to 
hf^i{t) = -0.996524 1 + 1.04747. 

unclear. The gap in the binary coding is smaller than in the polygonal description. 
This behaviour is generally expected by coarse graining the partition of Markov 
matrices, which in our case means grouping polygons in binary coded partition 
elements. The entropy of the Markov matrix in the binary encoding can be written 
in the following way 

Ht) = - E /"(^-) logr/, + (1 - r^^) log(l - r]^)] . (72) 

The entropy is bounded from above as 

h{t) < log(2) Aspiit , Aspiit = E H{r]^)fi{A^) , (73) 

ujeLt 



10" 



10" 



where Agput is the total area of elements of the partition Agput that are split by the 
action of the map 0. The number of split partition elements is N{t + 1) — N{t) and 
their total average area is Agput ~ {N{t + 1) — N{t))/N{t) ~ const t^^. Therefore 
we conclude that h{t) ~ const t^^ and this is numerically verified in Figure [26b . 
We calculate the matrix elements M^/^^^ at given t by a Monte-Carlo integration 
using iVjnc points uniformly sampled over the torus so that N{t)/N^c < 10~^. We 
empirically find that it suffices to have the relative error of the gap in observed time 
windows smaller than 5%. 

Next we consider a class of random Markov matrices {M G M"^" : J2j = 1} 
of dimension n ^ 1, which mimics known statistical properties of the Markov 
matrices obtained in the generic case of the triangle map system. We are interested 
in the asymptotics as n — oo. The individual matrix M in this class is a perturbed 
permutation matrix defined as 

M = R{id +D) , M e M![^" , (74) 
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where i? is a permutation matrix (jHeerstermaru . Il990l ) associated to a random inde- 

1 and D represent 



composable permutation (rj)"^o^5 with matrix elements Pij 
the perturbation in the expression fl74p given by 



Dij = USi,d^ - 5^j) , E [0, 1] 



(75) 



with {di)^!^ being a random permutation 
non-zero coefficients with an average ^ = 
between eigenvalues and the uni t circle in these class of Mar kov matrices using the 



There is only m = 0} <^ n 

— J2i 6- We estimate the average gap 



standard perturbation method (jCourant and Hilbertl . Il989l ). The indecomposible 



permutation matrix R is similar to down-shift permutation matrix 5* 



S = 0'^RO 



5, 



6^ 



i~j mod n,0 5 



(76) 



where O is just a product of transpositions. The matrix 5* can be diagonalized and 
its eigensystem is given by 

SF = FA, A = diag (e^^^)";^' , Fj,k = ■ (^7) 



Knowing that the eigensystem of the permutation matrix is RV = V A with V = 
OF = [vi]^ZQ, we find that amplitudes of eigenvalues Aj of matrix M in the ffist 
order perturbation are 

\Xi\ = \/l + 2mt{vjDvi} + \v}Dv,\^ , (78) 



and in the leading order depend only on the coefficients 

1 , 777 

|A,|^l--^e, = l--e. (79) 
n j n 



Following this formula (179|) the effective spectral gap for the Markov matrices of 
the triangle map should scale with time asymptotics as 

3 

A = l-max|Ai| , (80) 

« 2t 



where n = N{t) ~ Ct^, m = N{t + 1) — N{t) ~ 3Ct^ and, because the elements in 
the rows of matrix D can be equally large we take ^ = |- The results do not give 
the dependence of the gap in the polygonal and binary encoding of dynamics, but 
rather a correct time scaling of the distribution of the group of largest eigenvalues 
in both descriptions. 
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6 Stochastic models 



The deterministic model is found to be a hard system to analyse, and in particular 
to establish rigorous results. In order to avoid some of the most difficult technical 
problems, but still to obtain intuitively correct understanding of the dynamics of 
triangle maps we turn to stochastic models. The stochasticity is introduced in the 
position of the discontinuity of the triangle map (the cut), which is from analytical 
point of view, the most problematic. 



6.1 Random triangle map 

We introduce a stochastic version of the triangle map, called random triangle map, 
by randomising the position of the discontinuity and thereby preserving essential 
properties of the deterministic case e.g. the non-hyperbolic nature of the dynamics. 
The random triangle map is defined as the iteration formula 



where ^ [0, 1] is a u.d. random variable. In this system we again introduce the 
binary coding for a single realisation of the random sequence F = {E,t}tez* as it was 
introduced for the deterministic variant. 

We study the behaviour of partitions A(*)(r) with time t across realisations of the 
random sequence F. The number of partitions N{t, T) in binary coding at some t 
varies strongly with the realisation of the cuts F. Nevertheless we find a very simple 
behaviour of its mean N{t) and standard deviation cr^(t) in the generic case: 



N{t) := {Nit, F))r ~ ^^ := (iV(t, F)^)^ - N{tf ~ t'^ , (83) 



with Ci^2 being constants and (■)p denoting the average over all realisations of 
F. The numerically results are depicted in Figure [27b . where by fitting the data 
to the power law we find exponents p = 3.07 ± 0.02 and r = 1.84 ± 0.01. We 
may heuristically expect p = 3 as this is somehow the largest power we observe 
theoretically in sequences of such simple maps. We conclude that the number of 
elements of the partitions is asymptotically proportional to t^ similarly as in the 
deterministic variant. In addition the relative deviation a^/N = 0{t^~'^), where 
p > r, meaning that presented behaviour for the mean value is asymptotically 
stable. Furthermore we investigate the distribution of N{t, F) and we find a very 



Pt+i =Pt + a sgn(^t -qt)+ P mod 1 , 
qt+i = qt+ Pt+i mod 1 , 
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Fig. 27. The average number of the elements of the partitions N{t) as a function of time 
t (a) and distribution of the normaUsed number of elements of the partitions Ps{x,t) (b) 
in the generic case of random triangle map at q = 1/e and (3 = {V5-l)/2. 

good agreement with the Gaussian distribution for all t 



Ps{x,t) = {6{x- s(t,r)))i 



1 



exp 



for the rescaled number of elements of the partitions s(t, T) = {N(t, T)—N(t))/aNit)- 
This is demonstrated in Figure [27b. The form of the distribution and the same scal- 
ing of and with time t makes us to suspect that N{t, F) should be, for t ^ 1, 
approximately represented as a sum of ~ t independent stochastic variables. Next 
we take a look a the average distribution of areas of the elements of the partitions 
in the binary coding. We study the average distribution of areas of elements of the 
partitions 

n,{x,t) = / 5(x-MA))) , (85) 

\AGA(t){r) / p 



and average distribution of rescaled areas 

n^{x,t) = / 5(x-iV(t,F)/i(A))\ . (86) 

\AeA(*)(r) / p 

over the set of partitions. These two average distributions are not equivalent due to 
changes to the number elements of the partitions between different realisations of 
the cuts F. In generic case shown in Figure [28] we find that these distributions have 
similar functional dependence, which seems to have a simple asymptotics at large 
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Fig. 28. The average distribution of areas of the elements of partitions na,{x,t) (a) and of 
rescaled ares of the elements of partitions np{x, t) (b) in the generic case of triangle map 
at a = 1/e and /3 = (^/5 - l)/2. 



times reading 

lim t~^np{t^x,t) , ni,{x,t) ^ exp(-|0(x^)|) , (87) 

t — >oo 

The obtained dependence coincides with the dependence in the deterministic case 
using the binary coding, but due to the computational complexity of the problem 
we cannot establish the nature of the decay more precisely. 

Next we examine the mixing property in the stochastic model by studying the auto- 
correlation function of an observable f : T"^ M. with zero phase-space integral 
Jj2 dqdp f{q,p) = reading 

C{t,T)[f] = J dqdpfiq,p)fiqt{q,p,T),ptiq,p,T)), (88) 

T2 



where a trajectory with a given F and initial point {q,p) is denoted by 
{qt{q,p,T),pt{q,p,T)) for t G Z*. In particular, we focus on the first two central 
moments of the correlations w.r.t. the average over F i.e. C{t) := (C(t, F))p and 
cr^(t) := (C(t,F)^)p — C{t) . We performed numerical calculations of the two mo- 
ments using different parameters {a, (3). The results for the generic case of param- 
eters are shown in Figure [291 whereas for the non-generic case in Figure [301 In 
Figures [29k and [30k we see that the average autocorrelation C{t) in the generic and 
non-generic case decays exponentially 

C{t)=eM-\Om, (89) 
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Fig. 29. The average correlation C{t) (a) and its standard deviation (Tc{t) (b) for observ- 
ables f{q,p) = sin(27rg) + sin(27rp) (top) and f{q,p) = 9{q) (bottom) in random triangle 
model for a = 1/e and (5 = (\/5 — l)/2. The horizontal line in figures at the bottom 
represents the value of the statistical error of the method. 



and stretch-exponentially 
C(t)=exp(-|0(t^)|), 



(90) 



for large times t, respectively. The decay of average correlation is faster than alge- 
braic, which is expected in stochastic systems. In all tested cases we obtain in the 
random triangle map the same power law decay of the standard deviation of the 
autocorrelation reading 



(91) 



This is the same power-law decay as the one found for the correlation in generic case 
of the deterministic triangle map. This behaviour can be clearly seen in figures [29b 
and [30b. In Subsection 16.31 we clarify this behaviour by showing that the dispersion 
of correlation cr^(t) is inversely proportional to N{t), which in our case scales as 
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Fig. 30. The average correlation C(t) (a) and its standard deviation acit) (b) for ob- 
servables f{q,p) = 9{q) in random triangle model at a = 1/e and /3 = calculated by 
extrapolating results obtained by applying Simpson quadrature rule with = Ny = 500 
and averaging over 1600 realisations. 



6.2 Random triangle map acting on one- dimensional sets 



Here we would like to show that the model of random triangle map can be even 
simplified and still possesses some interesting features of the (random) triangle map. 

In the previous subsection we have analysed distributions of polygon/partition at- 
tributes within a stochastic model. The following random model indeed mimics the 
deterministic triangle map (P) using the polygonal description, if we assume random 
positions (x coordinates) of the cuts, relative to the positions of a given polygonal 
piece. In this subsection we show that such analysis can be performed analytically 
for a special case of a set of degenerate polygons, namely the line segments. We de- 
rive explicit results for the distribution of the length of line segments as a function 
of time. 

Let us take initial set to be a vertical line of certain length. Then after time t, the 
image is composed of many short lines, each of which has slope s = 1/t. So, the 
correct scaled slope variable is a constant Xs = st = 1 which is consistent with 
the scaling of the general case above fl37|) . Therefore we may continue investigating 
just one of the quantities, e.g. ■= h < 1, while the others are related to it, e.g. 
6y = h/t. Let P{h,t)dh be the probability that at time t a randomly chosen piece 
has a horizontal length between h — dh/2 and h + dh/2. Vertical size of this piece 
is h/t, so the total probability that this piece does not get cut in one iteration is 
1 — h. Here all pieces are taken with the statistical weight 1, so P{h,t) represents 
a distribution of the number of such pieces. Using elementary probabihty it is easy 
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to derive the following dynamic evolution integro-difference equation 

1 

P{h, t + 1) = (1 - h)P{h, t) + J dh' P{h', t) . (92) 



We are interested in a possible asymptotic scaling solutions of the above equation. 
Indeed, using the ansatz 

oo 

P{h, t) = ti){th) , j dt tp{t) = 1 , (93) 



for large t, the equation fl92l) is in the leading order in 1/t equivalent to an ordinary 
integro-differential equation 

oo 

^{x) + Xi^\x) + Xi^{x) = / dxtpix) , (94) 



which has a unique solution ip{x) = 6xp(— x) satisfying the appropriate boundary 
conditions, |'?/'(0)| < oo,ip{oo) = 0. Note that ip{x) is the scaled number density, 
■ip{x) = diVx.(x)/dX) so we have the cumulative distribution, 

iV,(x) = l-exp(-x), (95) 



The first remarkable result of this simplified analysis is the correct scaling variables 
Xy = f^Sy = t'^h, Xx = t^x = th. In the spirit of the presented simplified analysis 
we can also express the effective area of a polygon a to be proportional to 6y6x = 
XxXy/t^- Then the formula fl95l) immediately yields the cumulative distribution of 
the areas reading 

Naia) ^ 1 - exp(0(v^)) . (96) 



The second notable results is that the distributions ( l95i) and ( l96l) are essentially 
identical, or very close for the areas distribution, to the numerical results in the 
deterministic triangle map reported in Figure | 



6.3 Relation for the correlations decay 



In the following we present a heuristic derivation of the correlation decay in (ran- 
dom) triangle map using only its basic properties. Let us take two sets AjBtzT"^ 
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of positive measure and consider the correlations of characteristic functions of the 

sets, 

CaA^) = ^i{B n - ^l{B)M . (97) 



The total number of pieces to which A is cut up to time t, with average area (a) (t), 
see equations (IMI) . (158|1 and (183|) . is in all studied cases estimated by 



NAt) = ^44t ~ const for t > 1 

(a) (t) 



(98 



The number of pieces of A which map to a test set B, and hence contributes to the 
correlation function fl971). is 



Assuming that pieces (partitions) are uniformly and pseudo-randomly distributed 



over the set B, the expected standard deviation of Na,b is yN^^B, and an estimated 
value of the fluctuation of the correlation function reads 



^c^,e(i) = ^N^At) («) it) = ^/KA)fi{B) (a) (t) = Oit-'i) , (100) 

for t ^ 1, where we implicitly use that (C^,b(^)) is negligible in comparison to 
(Tc^ B (t) ■ This is precisely the scaling law which is observed in numerical explorations 
of the random triangle map. We believe that this elucidates the law of decay of 
correlations in the (deterministic) triangle map. 



7 Conclusions 



We presented both theoretical and numerical analysis of a two-parameter family 
of nonhyperbolic dynamical system related to the classical motion in a polygonally 
shaped billiard. We formulated several conjectures, and provided strong numerical 
evidence for their validity, concerning the ergodic behaviour of this family for differ- 
ent sets of parameters values. In particular, we believe that the system is strongly 
mixing (and thus ergodic) at least in a situation named generic, where both pa- 
rameters are mutually incommensurate irrationals so that in particular no periodic 
orbits exist. We also found a kind of mixing behaviour in the case where only the pa- 
rameter a controlling the strength of the cut is irrational whereas (3 is equal to zero. 
However, the two cases exhibit different mixing behaviours, diffusion properties and 
periodic orbit structure. 
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The dynamical properties are studied by means of two different symbolic encod- 
ing schemes based on propagating an initial set of polygons and following the two 
portions of phase space separated by the discontinuity lines, referred to as the polyg- 
onal and the binary description, respectively. Using these schemes we numerically 
deduced several asymptotic scaling properties - for some of which we provided a 
semi-heuristic analytical explanation. In particular we focused on the spectrum of 
the Markov transition matrix (in either of the two symbolic descriptions) near to 
the unit circle and on its spectral gap in relation to the correlation decay. The av- 
erage convergence of the group of largest eigenvalues to the unit circle as the size of 
the matrix increases is quantitatively explained by perturbative analysis of a simple 
random model. We have also briefly studied the Koopman operator associated to 
the map on a finite dimensional Hilbert space obtained by truncating the Fourier 
basis. There we found interesting scaling properties of the spectrum and fidelity 
with respect to the dimension of the space. 

Finally, we have introduced a stochastic version of the triangle map which allowed 
us to obtain some theoretical asymptotic scalings of the distribution of degenerate 
polygons, thus yielding the decay of correlations observed numerically in the 
generic case. 

This paper presents a bulk of numerical results, along with sparse analytical ar- 
guments, on the properties of the triangle map. From the mathematical point of 
view several important problems remain to be solved. Some of them, in particular 
the dynamical mechanism underlying the observed strongly mixing behaviour (for 
some parameter values) of this piecewise parabolic system is a challenging problem 
in crgodic theory whose solution may require the development of new ideas and 
techniques. 
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